rm(list=ls())
graphics.off()

install.packages("readxl")
install.packages("dplyr")
install.packages("tidyverse")
install.packages("ggpubr")
install.packages("rstatix")
install.packages("datarium")
install.packages("Biodem")
install.packages("multcomp")
install.packages("hrbrthemes")
install.packages("viridis")
install.packages("gapminder")
install.packages("devtools")
install.packages("stringr")
install.packages("patchwork")
install.packages("extrafont")
install.packages("gridExtra")
install.packages("gtable")
install.packages("ggtext")
install.packages("effsize")

library(ggplot2)
library(dplyr)
library(tidyverse)
library(ggpubr)
library(rstatix)
library(nlme)
library(matrixcalc)
library(Biodem)
library(stats)
library(multcomp)
library(hrbrthemes)
library(viridis)
library(gapminder)
library(devtools)
library(stringr)
library(patchwork)
library(gridExtra)
library(cowplot)
library(extrafont)
library(gtable)
library(ggtext)
library(effsize)

######################### START ###########################################

######### Read in data ##########

######################### START ###########################################

setwd("//nausers01/User/GOIKOG/Desktop/MatLab/Output/IntensityMeasures3")

data <- read.csv('//nausers01/User/GOIKOG/Desktop/Matlab/Output/IntensityMeasures3/NASATLX_DATA.csv', header=TRUE,sep=",")


#Tables
data_NASA_Test = bind_cols(data$Record.ID, data$Cognition.Easy.Total, data$Cognition.Challenging.Total, data$Cognition.Difficult.Total, data$Motor.Easy.Total, data$Motor.Challenging.Total, data$Motor.Difficult.Total)

data_NASA_ReTest = bind_cols(data$Record.ID, data$Cognition.Easy.2.Total, data$Cognition.Challenging.2.Total, data$Cognition.Difficult.2.Total, data$Motor.Easy.2.Total, data$Motor.Challenging.2.Total, data$Motor.Difficult.2.Total)

data_NASA_ReTest[10,] = data_NASA_Test[10,]
  
NASA = (data_NASA_Test+data_NASA_ReTest)/2


colnames(NASA) <- c("ID","aclow","bcchall","cchigh","dmlow","emchall","fmhigh")
colnames(data_NASA_Test) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")
colnames(data_NASA_ReTest) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")


###NASA 
NASA_long <- pivot_longer(NASA, cols=c("aclow","bcchall","cchigh","dmlow","emchall","fmhigh"))
div <- rep(c("Mental","Motor"), each = 3)
NASA_long$div <- rep(div,nrow(NASA)) 

NASA_long = split(NASA_long, NASA_long$div)

NASA_coglong = NASA_long$Mental[,1:3]

NASA_motlong =  NASA_long$Motor[,1:3]

###

plot <- ggplot(NASA_long, aes(name, value, ID, fill = name)) +
  geom_boxplot() +
  labs(x = 'Intensity Level',
       title = 'NASA',
       y = "NASA") +
  theme(plot.title = element_text(hjust = 0.5, size = 20),
        legend.title=element_blank(),
        panel.background = element_blank(),
        panel.grid.major.y = element_line(colour = "grey", size = 0.5),
        axis.title = element_text(size = 13),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        strip.text.x = element_text(size = 13)) +
  scale_x_discrete(label = element_blank()) +
  scale_y_continuous(limits = c(-0.015,0.223)) +
  guides(fill="none") +
  scale_fill_manual(values=c("steelblue4", "lightblue1", "cadetblue3", "steelblue4", "lightblue1", "cadetblue3")) +
  annotate(geom = "text", x = 0.83, y = -0.009, label = "Very Low", hjust = 0, vjust = 0, size = 3, color = "gray50") +  
  annotate(geom = "text", x = 0.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.009, label = "Challenging", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.85, y = -0.009, label = "Very High", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.87, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  facet_grid(. ~ div, scales = "free", space = "free")


#ggsave("HRVplot.png", plot, width = 10.5, height = 6.92)

###Statistics

###Test Normality
Norm.NASA.CogEasy = shapiro.test(NASA$aclow)
Norm.NASA.CogMiddle = shapiro.test(NASA$bcchall)
Norm.NASA.CogDifficult = shapiro.test(NASA$cchigh)

Norm.NASA.MotEasy = shapiro.test(NASA$dmlow)
Norm.NASA.MotMiddle = shapiro.test(NASA$emchall)
Norm.NASA.MotDifficult = shapiro.test(NASA$fmhigh)

###Perform Test and Post-Hoc
### Cognitive 
if (Norm.NASA.CogEasy$p.value > 0.05 & Norm.NASA.CogMiddle$p.value > 0.05 & Norm.NASA.CogDifficult$p.value > 0.05){
  print(anova_test(data=NASA_coglong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_coglong$value, NASA_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_coglong, value ~ name|ID))
  pairwise.wilcox.test(NASA_coglong$value, NASA_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}

###Motor
if (Norm.NASA.MotEasy$p.value > 0.05 & Norm.NASA.MotMiddle$p.value > 0.05 & Norm.NASA.MotDifficult$p.value > 0.05){
  print(anova_test(data=NASA_motlong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_motlong$value, NASA_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_motlong, value ~ name|ID))
  pairwise.wilcox.test(NASA_motlong$value, NASA_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}


#ICC analysis

## NASA

## Cognitive

### Easy
data_NASA_Cog_easy = data.table(data_NASA_Test$NASA_CognitionEasy, data_NASA_ReTest$NASA_CognitionEasy)
#data_NASA_easy_RRn = data_NASA_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Cog_easy =ICC(data_NASA_Cog_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Cog_easy<- bland.altman.plot(data_NASA_Test$NASA_CognitionEasy,data_NASA_ReTest$NASA_CognitionEasy, main="Test-Retest NASA Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Cog_easy)

###Middle
data_NASA_Cog_middle = data.table(data_NASA_Test$NASA_CognitionOptimal, data_NASA_ReTest$NASA_CognitionOptimal)
#data_NASA_Cog_middle = data_NASA_Cog_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Cog_middle =ICC(data_NASA_Cog_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Cog_middle <-bland.altman.plot(data_NASA_Test$NASA_CognitionOptimal,data_NASA_ReTest$NASA_CognitionOptimal, main="Test-Retest NASA Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Cog_middle)


###Difficult
data_NASA_Cog_difficult = data.table(data_NASA_Test$NASA_CognitionDifficult, data_NASA_ReTest$NASA_CognitionDifficult)
#data_NASA_Cog_difficult = data_NASA_Cog_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Cog_difficult =ICC(data_NASA_Cog_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Cog_difficult<- bland.altman.plot(data_NASA_Test$NASA_CognitionDifficult,data_NASA_ReTest$NASA_CognitionDifficult, main="Test-Retest NASA Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Cog_difficult)

## Motor

### Easy
data_NASA_Mot_easy = data.table(data_NASA_Test$NASA_MotorEasy, data_NASA_ReTest$NASA_MotorEasy)
#data_NASA_easy_RRn = data_NASA_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mot_easy =ICC(data_NASA_Mot_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mot_easy<- bland.altman.plot(data_NASA_Test$NASA_MotorEasy,data_NASA_ReTest$NASA_MotorEasy, main="Test-Retest NASA Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Mot_easy)

###Middle
data_NASA_Mot_middle = data.table(data_NASA_Test$NASA_MotorOptimal, data_NASA_ReTest$NASA_MotorOptimal)
#data_NASA_Mot_middle = data_NASA_Mot_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mot_middle =ICC(data_NASA_Mot_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mot_middle <-bland.altman.plot(data_NASA_Test$NASA_MotorOptimal,data_NASA_ReTest$NASA_MotorOptimal, main="Test-Retest NASA Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Mot_middle)


###Difficult
data_NASA_Mot_difficult = data.table(data_NASA_Test$NASA_MotorDifficult, data_NASA_ReTest$NASA_MotorDifficult)
#data_NASA_Mot_difficult = data_NASA_Mot_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mot_difficult =ICC(data_NASA_Mot_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mot_difficult<- bland.altman.plot(data_NASA_Test$NASA_MotorDifficult,data_NASA_ReTest$NASA_MotorDifficult, main="Test-Retest NASA Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Mot_difficult)


########################################################################################################################################################
########################################################################################################################################################
########################################################################################################################################################



### NASA by DIMENSION


###NASA_Mental

data_NASA_Mental_Test = bind_cols(data$Record.ID, data$Mental.Demand, data$Mental.Demand.1 ,data$Mental.Demand.2 ,data$Mental.Demand.3 ,data$Mental.Demand.4 , data$Mental.Demand.5)

data_NASA_Mental_ReTest = bind_cols(data$Record.ID, data$Mental.Demand.6, data$Mental.Demand.7 ,data$Mental.Demand.8 ,data$Mental.Demand.9 ,data$Mental.Demand.10 , data$Mental.Demand.11)

data_NASA_Mental_ReTest[10,] = data_NASA_Mental_Test[10,]

NASA_Mental = (data_NASA_Mental_Test+data_NASA_Mental_ReTest)/2

colnames(NASA_Mental) <- c("ID","aclow","bcchall","cchigh","dmlow","emchall","fmhigh")
colnames(data_NASA_Mental_Test) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")
colnames(data_NASA_Mental_ReTest) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")



NASA_Mental_long <- pivot_longer(NASA_Mental, cols=c("aclow","bcchall","cchigh","dmlow","emchall","fmhigh"))
div <- rep(c("Mental","Motor"), each = 3)
NASA_Mental_long$div <- rep(div,nrow(NASA)) 

NASA_Mental_long = split(NASA_Mental_long, NASA_Mental_long$div)

NASA_Mental_coglong = NASA_Mental_long$Mental[,1:3]

NASA_Mental_motlong =  NASA_Mental_long$Motor[,1:3]

###

plot <- ggplot(NASA_Mental_long, aes(name, value, ID, fill = name)) +
  geom_boxplot() +
  labs(x = 'Intensity Level',
       title = 'NASA_Mental',
       y = "NASA_Mental") +
  theme(plot.title = element_text(hjust = 0.5, size = 20),
        legend.title=element_blank(),
        panel.background = element_blank(),
        panel.grid.major.y = element_line(colour = "grey", size = 0.5),
        axis.title = element_text(size = 13),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        strip.text.x = element_text(size = 13)) +
  scale_x_discrete(label = element_blank()) +
  scale_y_continuous(limits = c(-0.015,0.223)) +
  guides(fill="none") +
  scale_fill_manual(values=c("steelblue4", "lightblue1", "cadetblue3", "steelblue4", "lightblue1", "cadetblue3")) +
  annotate(geom = "text", x = 0.83, y = -0.009, label = "Very Low", hjust = 0, vjust = 0, size = 3, color = "gray50") +  
  annotate(geom = "text", x = 0.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.009, label = "Challenging", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.85, y = -0.009, label = "Very High", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.87, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  facet_grid(. ~ div, scales = "free", space = "free")


#ggsave("HRVplot.png", plot, width = 10.5, height = 6.92)

###Statistics

###Test Normality
Norm.NASA_Mental.CogEasy = shapiro.test(NASA_Mental$aclow)
Norm.NASA_Mental.CogMiddle = shapiro.test(NASA_Mental$bcchall)
Norm.NASA_Mental.CogDifficult = shapiro.test(NASA_Mental$cchigh)

Norm.NASA_Mental.MotEasy = shapiro.test(NASA_Mental$dmlow)
Norm.NASA_Mental.MotMiddle = shapiro.test(NASA_Mental$emchall)
Norm.NASA_Mental.MotDifficult = shapiro.test(NASA_Mental$fmhigh)

###Perform Test and Post-Hoc
### Cognitive 
if (Norm.NASA_Mental.CogEasy$p.value > 0.05 & Norm.NASA_Mental.CogMiddle$p.value > 0.05 & Norm.NASA_Mental.CogDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Mental_coglong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Mental_coglong$value, NASA_Mental_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Mental_coglong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Mental_coglong$value, NASA_Mental_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}

###Motor
if (Norm.NASA_Mental.MotEasy$p.value > 0.05 & Norm.NASA_Mental.MotMiddle$p.value > 0.05 & Norm.NASA_Mental.MotDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Mental_motlong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Mental_motlong$value, NASA_Mental_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Mental_motlong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Mental_motlong$value, NASA_Mental_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}


#ICC analysis

## NASA_Mental

## Cognitive

### Easy
data_NASA_Mental_Cog_easy = data.table(data_NASA_Mental_Test$NASA_CognitionEasy, data_NASA_Mental_ReTest$NASA_CognitionEasy)
#data_NASA_Mental_easy_RRn = data_NASA_Mental_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mental_Cog_easy =ICC(data_NASA_Mental_Cog_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mental_Cog_easy<- bland.altman.plot(data_NASA_Mental_Test$NASA_CognitionEasy,data_NASA_Mental_ReTest$NASA_CognitionEasy, main="Test-Retest NASA_Mental Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Mental_Cog_easy)

###Middle
data_NASA_Mental_Cog_middle = data.table(data_NASA_Mental_Test$NASA_CognitionOptimal, data_NASA_Mental_ReTest$NASA_CognitionOptimal)
#data_NASA_Mental_Cog_middle = data_NASA_Mental_Cog_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mental_Cog_middle =ICC(data_NASA_Mental_Cog_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mental_Cog_middle <-bland.altman.plot(data_NASA_Mental_Test$NASA_CognitionOptimal,data_NASA_Mental_ReTest$NASA_CognitionOptimal, main="Test-Retest NASA_Mental Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Mental_Cog_middle)


###Difficult
data_NASA_Mental_Cog_difficult = data.table(data_NASA_Mental_Test$NASA_CognitionDifficult, data_NASA_Mental_ReTest$NASA_CognitionDifficult)
#data_NASA_Mental_Cog_difficult = data_NASA_Mental_Cog_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mental_Cog_difficult =ICC(data_NASA_Mental_Cog_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mental_Cog_difficult<- bland.altman.plot(data_NASA_Mental_Test$NASA_CognitionDifficult,data_NASA_Mental_ReTest$NASA_CognitionDifficult, main="Test-Retest NASA_Mental Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Mental_Cog_difficult)

## Motor

### Easy
data_NASA_Mental_Mot_easy = data.table(data_NASA_Mental_Test$NASA_MotorEasy, data_NASA_Mental_ReTest$NASA_MotorEasy)
#data_NASA_Mental_easy_RRn = data_NASA_Mental_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mental_Mot_easy =ICC(data_NASA_Mental_Mot_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mental_Mot_easy<- bland.altman.plot(data_NASA_Mental_Test$NASA_MotorEasy,data_NASA_Mental_ReTest$NASA_MotorEasy, main="Test-Retest NASA_Mental Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Mental_Mot_easy)

###Middle
data_NASA_Mental_Mot_middle = data.table(data_NASA_Mental_Test$NASA_MotorOptimal, data_NASA_Mental_ReTest$NASA_MotorOptimal)
#data_NASA_Mental_Mot_middle = data_NASA_Mental_Mot_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mental_Mot_middle =ICC(data_NASA_Mental_Mot_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mental_Mot_middle <-bland.altman.plot(data_NASA_Mental_Test$NASA_MotorOptimal,data_NASA_Mental_ReTest$NASA_MotorOptimal, main="Test-Retest NASA_Mental Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Mental_Mot_middle)


###Difficult
data_NASA_Mental_Mot_difficult = data.table(data_NASA_Mental_Test$NASA_MotorDifficult, data_NASA_Mental_ReTest$NASA_MotorDifficult)
#data_NASA_Mental_Mot_difficult = data_NASA_Mental_Mot_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Mental_Mot_difficult =ICC(data_NASA_Mental_Mot_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Mental_Mot_difficult<- bland.altman.plot(data_NASA_Mental_Test$NASA_MotorDifficult,data_NASA_Mental_ReTest$NASA_MotorDifficult, main="Test-Retest NASA_Mental Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Mental_Mot_difficult)




########################################################################################################################################################


###NASA_Temporal

data_NASA_Temporal_Test = bind_cols(data$Record.ID, data$Temporal.Demand, data$Temporal.Demand.1 ,data$Temporal.Demand.2 ,data$Temporal.Demand.3 ,data$Temporal.Demand.4 , data$Temporal.Demand.5)

data_NASA_Temporal_ReTest = bind_cols(data$Record.ID, data$Temporal.Demand.6, data$Temporal.Demand.7 ,data$Temporal.Demand.8 ,data$Temporal.Demand.9 ,data$Temporal.Demand.10 , data$Temporal.Demand.11)

data_NASA_Temporal_ReTest[10,] = data_NASA_Temporal_Test[10,]

NASA_Temporal = (data_NASA_Temporal_Test+data_NASA_Temporal_ReTest)/2

colnames(NASA_Temporal) <- c("ID","aclow","bcchall","cchigh","dmlow","emchall","fmhigh")
colnames(data_NASA_Temporal_Test) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")
colnames(data_NASA_Temporal_ReTest) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")



NASA_Temporal_long <- pivot_longer(NASA_Temporal, cols=c("aclow","bcchall","cchigh","dmlow","emchall","fmhigh"))
div <- rep(c("Mental","Motor"), each = 3)
NASA_Temporal_long$div <- rep(div,nrow(NASA_Temporal)) 

NASA_Temporal_long = split(NASA_Temporal_long, NASA_Temporal_long$div)

NASA_Temporal_coglong = NASA_Temporal_long$Mental[,1:3]

NASA_Temporal_motlong =  NASA_Temporal_long$Motor[,1:3]

###

plot <- ggplot(NASA_Temporal_long, aes(name, value, ID, fill = name)) +
  geom_boxplot() +
  labs(x = 'Intensity Level',
       title = 'NASA_Temporal',
       y = "NASA_Temporal") +
  theme(plot.title = element_text(hjust = 0.5, size = 20),
        legend.title=element_blank(),
        panel.background = element_blank(),
        panel.grid.major.y = element_line(colour = "grey", size = 0.5),
        axis.title = element_text(size = 13),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        strip.text.x = element_text(size = 13)) +
  scale_x_discrete(label = element_blank()) +
  scale_y_continuous(limits = c(-0.015,0.223)) +
  guides(fill="none") +
  scale_fill_manual(values=c("steelblue4", "lightblue1", "cadetblue3", "steelblue4", "lightblue1", "cadetblue3")) +
  annotate(geom = "text", x = 0.83, y = -0.009, label = "Very Low", hjust = 0, vjust = 0, size = 3, color = "gray50") +  
  annotate(geom = "text", x = 0.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.009, label = "Challenging", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.85, y = -0.009, label = "Very High", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.87, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  facet_grid(. ~ div, scales = "free", space = "free")


#ggsave("HRVplot.png", plot, width = 10.5, height = 6.92)

###Statistics

###Test Normality
Norm.NASA_Temporal.CogEasy = shapiro.test(NASA_Temporal$aclow)
Norm.NASA_Temporal.CogMiddle = shapiro.test(NASA_Temporal$bcchall)
Norm.NASA_Temporal.CogDifficult = shapiro.test(NASA_Temporal$cchigh)

Norm.NASA_Temporal.MotEasy = shapiro.test(NASA_Temporal$dmlow)
Norm.NASA_Temporal.MotMiddle = shapiro.test(NASA_Temporal$emchall)
Norm.NASA_Temporal.MotDifficult = shapiro.test(NASA_Temporal$fmhigh)

###Perform Test and Post-Hoc
### Cognitive 
if (Norm.NASA_Temporal.CogEasy$p.value > 0.05 & Norm.NASA_Temporal.CogMiddle$p.value > 0.05 & Norm.NASA_Temporal.CogDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Temporal_coglong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Temporal_coglong$value, NASA_Temporal_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Temporal_coglong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Temporal_coglong$value, NASA_Temporal_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}

###Motor
if (Norm.NASA_Temporal.MotEasy$p.value > 0.05 & Norm.NASA_Temporal.MotMiddle$p.value > 0.05 & Norm.NASA_Temporal.MotDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Temporal_motlong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Temporal_motlong$value, NASA_Temporal_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Temporal_motlong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Temporal_motlong$value, NASA_Temporal_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}


#ICC analysis

## NASA_Temporal

## Cognitive

### Easy
data_NASA_Temporal_Cog_easy = data.table(data_NASA_Temporal_Test$NASA_CognitionEasy, data_NASA_Temporal_ReTest$NASA_CognitionEasy)
#data_NASA_Temporal_easy_RRn = data_NASA_Temporal_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Temporal_Cog_easy =ICC(data_NASA_Temporal_Cog_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Temporal_Cog_easy<- bland.altman.plot(data_NASA_Temporal_Test$NASA_CognitionEasy,data_NASA_Temporal_ReTest$NASA_CognitionEasy, main="Test-Retest NASA_Temporal Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Temporal_Cog_easy)

###Middle
data_NASA_Temporal_Cog_middle = data.table(data_NASA_Temporal_Test$NASA_CognitionOptimal, data_NASA_Temporal_ReTest$NASA_CognitionOptimal)
#data_NASA_Temporal_Cog_middle = data_NASA_Temporal_Cog_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Temporal_Cog_middle =ICC(data_NASA_Temporal_Cog_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Temporal_Cog_middle <-bland.altman.plot(data_NASA_Temporal_Test$NASA_CognitionOptimal,data_NASA_Temporal_ReTest$NASA_CognitionOptimal, main="Test-Retest NASA_Temporal Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Temporal_Cog_middle)


###Difficult
data_NASA_Temporal_Cog_difficult = data.table(data_NASA_Temporal_Test$NASA_CognitionDifficult, data_NASA_Temporal_ReTest$NASA_CognitionDifficult)
#data_NASA_Temporal_Cog_difficult = data_NASA_Temporal_Cog_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Temporal_Cog_difficult =ICC(data_NASA_Temporal_Cog_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Temporal_Cog_difficult<- bland.altman.plot(data_NASA_Temporal_Test$NASA_CognitionDifficult,data_NASA_Temporal_ReTest$NASA_CognitionDifficult, main="Test-Retest NASA_Temporal Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Temporal_Cog_difficult)

## Motor

### Easy
data_NASA_Temporal_Mot_easy = data.table(data_NASA_Temporal_Test$NASA_MotorEasy, data_NASA_Temporal_ReTest$NASA_MotorEasy)
#data_NASA_Temporal_easy_RRn = data_NASA_Temporal_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Temporal_Mot_easy =ICC(data_NASA_Temporal_Mot_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Temporal_Mot_easy<- bland.altman.plot(data_NASA_Temporal_Test$NASA_MotorEasy,data_NASA_Temporal_ReTest$NASA_MotorEasy, main="Test-Retest NASA_Temporal Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Temporal_Mot_easy)

###Middle
data_NASA_Temporal_Mot_middle = data.table(data_NASA_Temporal_Test$NASA_MotorOptimal, data_NASA_Temporal_ReTest$NASA_MotorOptimal)
#data_NASA_Temporal_Mot_middle = data_NASA_Temporal_Mot_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Temporal_Mot_middle =ICC(data_NASA_Temporal_Mot_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Temporal_Mot_middle <-bland.altman.plot(data_NASA_Temporal_Test$NASA_MotorOptimal,data_NASA_Temporal_ReTest$NASA_MotorOptimal, main="Test-Retest NASA_Temporal Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Temporal_Mot_middle)


###Difficult
data_NASA_Temporal_Mot_difficult = data.table(data_NASA_Temporal_Test$NASA_MotorDifficult, data_NASA_Temporal_ReTest$NASA_MotorDifficult)
#data_NASA_Temporal_Mot_difficult = data_NASA_Temporal_Mot_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Temporal_Mot_difficult =ICC(data_NASA_Temporal_Mot_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Temporal_Mot_difficult<- bland.altman.plot(data_NASA_Temporal_Test$NASA_MotorDifficult,data_NASA_Temporal_ReTest$NASA_MotorDifficult, main="Test-Retest NASA_Temporal Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Temporal_Mot_difficult)




########################################################################################################################################################


###NASA_Physical

data_NASA_Physical_Test = bind_cols(data$Record.ID, data$Physical.Demand, data$Physical.Demand.1 ,data$Physical.Demand.2 ,data$Physical.Demand.3 ,data$Physical.Demand.4 , data$Physical.Demand.5)

data_NASA_Physical_ReTest = bind_cols(data$Record.ID, data$Physical.Demand.6, data$Physical.Demand.7 ,data$Physical.Demand.8 ,data$Physical.Demand.9 ,data$Physical.Demand.10 , data$Physical.Demand.11)

data_NASA_Physical_ReTest[10,] = data_NASA_Physical_Test[10,]

NASA_Physical = (data_NASA_Physical_Test+data_NASA_Physical_ReTest)/2

colnames(NASA_Physical) <- c("ID","aclow","bcchall","cchigh","dmlow","emchall","fmhigh")
colnames(data_NASA_Physical_Test) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")
colnames(data_NASA_Physical_ReTest) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")



NASA_Physical_long <- pivot_longer(NASA_Physical, cols=c("aclow","bcchall","cchigh","dmlow","emchall","fmhigh"))
div <- rep(c("Mental","Motor"), each = 3)
NASA_Physical_long$div <- rep(div,nrow(NASA_Physical)) 

NASA_Physical_long = split(NASA_Physical_long, NASA_Physical_long$div)

NASA_Physical_coglong = NASA_Physical_long$Mental[,1:3]

NASA_Physical_motlong =  NASA_Physical_long$Motor[,1:3]

###

plot <- ggplot(NASA_Physical_long, aes(name, value, ID, fill = name)) +
  geom_boxplot() +
  labs(x = 'Intensity Level',
       title = 'NASA_Physical',
       y = "NASA_Physical") +
  theme(plot.title = element_text(hjust = 0.5, size = 20),
        legend.title=element_blank(),
        panel.background = element_blank(),
        panel.grid.major.y = element_line(colour = "grey", size = 0.5),
        axis.title = element_text(size = 13),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        strip.text.x = element_text(size = 13)) +
  scale_x_discrete(label = element_blank()) +
  scale_y_continuous(limits = c(-0.015,0.223)) +
  guides(fill="none") +
  scale_fill_manual(values=c("steelblue4", "lightblue1", "cadetblue3", "steelblue4", "lightblue1", "cadetblue3")) +
  annotate(geom = "text", x = 0.83, y = -0.009, label = "Very Low", hjust = 0, vjust = 0, size = 3, color = "gray50") +  
  annotate(geom = "text", x = 0.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.009, label = "Challenging", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.85, y = -0.009, label = "Very High", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.87, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  facet_grid(. ~ div, scales = "free", space = "free")


#ggsave("HRVplot.png", plot, width = 10.5, height = 6.92)

###Statistics

###Test Normality
Norm.NASA_Physical.CogEasy = shapiro.test(NASA_Physical$aclow)
Norm.NASA_Physical.CogMiddle = shapiro.test(NASA_Physical$bcchall)
Norm.NASA_Physical.CogDifficult = shapiro.test(NASA_Physical$cchigh)

Norm.NASA_Physical.MotEasy = shapiro.test(NASA_Physical$dmlow)
Norm.NASA_Physical.MotMiddle = shapiro.test(NASA_Physical$emchall)
Norm.NASA_Physical.MotDifficult = shapiro.test(NASA_Physical$fmhigh)

###Perform Test and Post-Hoc
### Cognitive 
if (Norm.NASA_Physical.CogEasy$p.value > 0.05 & Norm.NASA_Physical.CogMiddle$p.value > 0.05 & Norm.NASA_Physical.CogDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Physical_coglong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Physical_coglong$value, NASA_Physical_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Physical_coglong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Physical_coglong$value, NASA_Physical_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}

###Motor
if (Norm.NASA_Physical.MotEasy$p.value > 0.05 & Norm.NASA_Physical.MotMiddle$p.value > 0.05 & Norm.NASA_Physical.MotDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Physical_motlong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Physical_motlong$value, NASA_Physical_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Physical_motlong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Physical_motlong$value, NASA_Physical_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}


#ICC analysis

## NASA_Physical

## Cognitive

### Easy
data_NASA_Physical_Cog_easy = data.table(data_NASA_Physical_Test$NASA_CognitionEasy, data_NASA_Physical_ReTest$NASA_CognitionEasy)
#data_NASA_Physical_easy_RRn = data_NASA_Physical_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Physical_Cog_easy =ICC(data_NASA_Physical_Cog_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Physical_Cog_easy<- bland.altman.plot(data_NASA_Physical_Test$NASA_CognitionEasy,data_NASA_Physical_ReTest$NASA_CognitionEasy, main="Test-Retest NASA_Physical Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Physical_Cog_easy)

###Middle
data_NASA_Physical_Cog_middle = data.table(data_NASA_Physical_Test$NASA_CognitionOptimal, data_NASA_Physical_ReTest$NASA_CognitionOptimal)
#data_NASA_Physical_Cog_middle = data_NASA_Physical_Cog_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Physical_Cog_middle =ICC(data_NASA_Physical_Cog_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Physical_Cog_middle <-bland.altman.plot(data_NASA_Physical_Test$NASA_CognitionOptimal,data_NASA_Physical_ReTest$NASA_CognitionOptimal, main="Test-Retest NASA_Physical Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Physical_Cog_middle)


###Difficult
data_NASA_Physical_Cog_difficult = data.table(data_NASA_Physical_Test$NASA_CognitionDifficult, data_NASA_Physical_ReTest$NASA_CognitionDifficult)
#data_NASA_Physical_Cog_difficult = data_NASA_Physical_Cog_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Physical_Cog_difficult =ICC(data_NASA_Physical_Cog_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Physical_Cog_difficult<- bland.altman.plot(data_NASA_Physical_Test$NASA_CognitionDifficult,data_NASA_Physical_ReTest$NASA_CognitionDifficult, main="Test-Retest NASA_Physical Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Physical_Cog_difficult)

## Motor

### Easy
data_NASA_Physical_Mot_easy = data.table(data_NASA_Physical_Test$NASA_MotorEasy, data_NASA_Physical_ReTest$NASA_MotorEasy)
#data_NASA_Physical_easy_RRn = data_NASA_Physical_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Physical_Mot_easy =ICC(data_NASA_Physical_Mot_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Physical_Mot_easy<- bland.altman.plot(data_NASA_Physical_Test$NASA_MotorEasy,data_NASA_Physical_ReTest$NASA_MotorEasy, main="Test-Retest NASA_Physical Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Physical_Mot_easy)

###Middle
data_NASA_Physical_Mot_middle = data.table(data_NASA_Physical_Test$NASA_MotorOptimal, data_NASA_Physical_ReTest$NASA_MotorOptimal)
#data_NASA_Physical_Mot_middle = data_NASA_Physical_Mot_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Physical_Mot_middle =ICC(data_NASA_Physical_Mot_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Physical_Mot_middle <-bland.altman.plot(data_NASA_Physical_Test$NASA_MotorOptimal,data_NASA_Physical_ReTest$NASA_MotorOptimal, main="Test-Retest NASA_Physical Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Physical_Mot_middle)


###Difficult
data_NASA_Physical_Mot_difficult = data.table(data_NASA_Physical_Test$NASA_MotorDifficult, data_NASA_Physical_ReTest$NASA_MotorDifficult)
#data_NASA_Physical_Mot_difficult = data_NASA_Physical_Mot_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Physical_Mot_difficult =ICC(data_NASA_Physical_Mot_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Physical_Mot_difficult<- bland.altman.plot(data_NASA_Physical_Test$NASA_MotorDifficult,data_NASA_Physical_ReTest$NASA_MotorDifficult, main="Test-Retest NASA_Physical Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Physical_Mot_difficult)




########################################################################################################################################################


###NASA_Performance

data_NASA_Performance_Test = bind_cols(data$Record.ID, data$Performance, data$Performance.1 ,data$Performance.2 ,data$Performance.3 ,data$Performance.4 , data$Performance.5)

data_NASA_Performance_ReTest = bind_cols(data$Record.ID, data$Performance.6, data$Performance.7 ,data$Performance.8 ,data$Performance.9 ,data$Performance.10 , data$Performance.11)

data_NASA_Performance_ReTest[10,] = data_NASA_Performance_Test[10,]

NASA_Performance = (data_NASA_Performance_Test+data_NASA_Performance_ReTest)/2

colnames(NASA_Performance) <- c("ID","aclow","bcchall","cchigh","dmlow","emchall","fmhigh")
colnames(data_NASA_Performance_Test) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")
colnames(data_NASA_Performance_ReTest) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")



NASA_Performance_long <- pivot_longer(NASA_Performance, cols=c("aclow","bcchall","cchigh","dmlow","emchall","fmhigh"))
div <- rep(c("Mental","Motor"), each = 3)
NASA_Performance_long$div <- rep(div,nrow(NASA_Performance)) 

NASA_Performance_long = split(NASA_Performance_long, NASA_Performance_long$div)

NASA_Performance_coglong = NASA_Performance_long$Mental[,1:3]

NASA_Performance_motlong =  NASA_Performance_long$Motor[,1:3]

###

plot <- ggplot(NASA_Performance_long, aes(name, value, ID, fill = name)) +
  geom_boxplot() +
  labs(x = 'Intensity Level',
       title = 'NASA_Performance',
       y = "NASA_Performance") +
  theme(plot.title = element_text(hjust = 0.5, size = 20),
        legend.title=element_blank(),
        panel.background = element_blank(),
        panel.grid.major.y = element_line(colour = "grey", size = 0.5),
        axis.title = element_text(size = 13),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        strip.text.x = element_text(size = 13)) +
  scale_x_discrete(label = element_blank()) +
  scale_y_continuous(limits = c(-0.015,0.223)) +
  guides(fill="none") +
  scale_fill_manual(values=c("steelblue4", "lightblue1", "cadetblue3", "steelblue4", "lightblue1", "cadetblue3")) +
  annotate(geom = "text", x = 0.83, y = -0.009, label = "Very Low", hjust = 0, vjust = 0, size = 3, color = "gray50") +  
  annotate(geom = "text", x = 0.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.009, label = "Challenging", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.85, y = -0.009, label = "Very High", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.87, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  facet_grid(. ~ div, scales = "free", space = "free")


#ggsave("HRVplot.png", plot, width = 10.5, height = 6.92)

###Statistics

###Test Normality
Norm.NASA_Performance.CogEasy = shapiro.test(NASA_Performance$aclow)
Norm.NASA_Performance.CogMiddle = shapiro.test(NASA_Performance$bcchall)
Norm.NASA_Performance.CogDifficult = shapiro.test(NASA_Performance$cchigh)

Norm.NASA_Performance.MotEasy = shapiro.test(NASA_Performance$dmlow)
Norm.NASA_Performance.MotMiddle = shapiro.test(NASA_Performance$emchall)
Norm.NASA_Performance.MotDifficult = shapiro.test(NASA_Performance$fmhigh)

###Perform Test and Post-Hoc
### Cognitive 
if (Norm.NASA_Performance.CogEasy$p.value > 0.05 & Norm.NASA_Performance.CogMiddle$p.value > 0.05 & Norm.NASA_Performance.CogDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Performance_coglong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Performance_coglong$value, NASA_Performance_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Performance_coglong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Performance_coglong$value, NASA_Performance_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}

###Motor
### Cognitive 
if (Norm.NASA_Performance.MotEasy$p.value > 0.05 & Norm.NASA_Performance.MotMiddle$p.value > 0.05 & Norm.NASA_Performance.MotDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Performance_motlong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Performance_motlong$value, NASA_Performance_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Performance_motlong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Performance_motlong$value, NASA_Performance_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}


#ICC analysis

## NASA_Performance

## Cognitive

### Easy
data_NASA_Performance_Cog_easy = data.table(data_NASA_Performance_Test$NASA_CognitionEasy, data_NASA_Performance_ReTest$NASA_CognitionEasy)
#data_NASA_Performance_easy_RRn = data_NASA_Performance_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Performance_Cog_easy =ICC(data_NASA_Performance_Cog_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Performance_Cog_easy<- bland.altman.plot(data_NASA_Performance_Test$NASA_CognitionEasy,data_NASA_Performance_ReTest$NASA_CognitionEasy, main="Test-Retest NASA_Performance Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Performance_Cog_easy)

###Middle
data_NASA_Performance_Cog_middle = data.table(data_NASA_Performance_Test$NASA_CognitionOptimal, data_NASA_Performance_ReTest$NASA_CognitionOptimal)
#data_NASA_Performance_Cog_middle = data_NASA_Performance_Cog_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Performance_Cog_middle =ICC(data_NASA_Performance_Cog_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Performance_Cog_middle <-bland.altman.plot(data_NASA_Performance_Test$NASA_CognitionOptimal,data_NASA_Performance_ReTest$NASA_CognitionOptimal, main="Test-Retest NASA_Performance Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Performance_Cog_middle)


###Difficult
data_NASA_Performance_Cog_difficult = data.table(data_NASA_Performance_Test$NASA_CognitionDifficult, data_NASA_Performance_ReTest$NASA_CognitionDifficult)
#data_NASA_Performance_Cog_difficult = data_NASA_Performance_Cog_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Performance_Cog_difficult =ICC(data_NASA_Performance_Cog_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Performance_Cog_difficult<- bland.altman.plot(data_NASA_Performance_Test$NASA_CognitionDifficult,data_NASA_Performance_ReTest$NASA_CognitionDifficult, main="Test-Retest NASA_Performance Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Performance_Cog_difficult)

## Motor

### Easy
data_NASA_Performance_Mot_easy = data.table(data_NASA_Performance_Test$NASA_MotorEasy, data_NASA_Performance_ReTest$NASA_MotorEasy)
#data_NASA_Performance_easy_RRn = data_NASA_Performance_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Performance_Mot_easy =ICC(data_NASA_Performance_Mot_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Performance_Mot_easy<- bland.altman.plot(data_NASA_Performance_Test$NASA_MotorEasy,data_NASA_Performance_ReTest$NASA_MotorEasy, main="Test-Retest NASA_Performance Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Performance_Mot_easy)

###Middle
data_NASA_Performance_Mot_middle = data.table(data_NASA_Performance_Test$NASA_MotorOptimal, data_NASA_Performance_ReTest$NASA_MotorOptimal)
#data_NASA_Performance_Mot_middle = data_NASA_Performance_Mot_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Performance_Mot_middle =ICC(data_NASA_Performance_Mot_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Performance_Mot_middle <-bland.altman.plot(data_NASA_Performance_Test$NASA_MotorOptimal,data_NASA_Performance_ReTest$NASA_MotorOptimal, main="Test-Retest NASA_Performance Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Performance_Mot_middle)


###Difficult
data_NASA_Performance_Mot_difficult = data.table(data_NASA_Performance_Test$NASA_MotorDifficult, data_NASA_Performance_ReTest$NASA_MotorDifficult)
#data_NASA_Performance_Mot_difficult = data_NASA_Performance_Mot_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Performance_Mot_difficult =ICC(data_NASA_Performance_Mot_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Performance_Mot_difficult<- bland.altman.plot(data_NASA_Performance_Test$NASA_MotorDifficult,data_NASA_Performance_ReTest$NASA_MotorDifficult, main="Test-Retest NASA_Performance Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Performance_Mot_difficult)





########################################################################################################################################################


###NASA_Effort

data_NASA_Effort_Test = bind_cols(data$Record.ID, data$Effort, data$Effort.1 ,data$Effort.2 ,data$Effort.3 ,data$Effort.4 , data$Effort.5)

data_NASA_Effort_ReTest = bind_cols(data$Record.ID, data$Effort.6, data$Effort.7 ,data$Effort.8 ,data$Effort.9 ,data$Effort.10 , data$Effort.11)

data_NASA_Effort_ReTest[10,] = data_NASA_Effort_Test[10,]

NASA_Effort = (data_NASA_Effort_Test+data_NASA_Effort_ReTest)/2

colnames(NASA_Effort) <- c("ID","aclow","bcchall","cchigh","dmlow","emchall","fmhigh")
colnames(data_NASA_Effort_Test) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")
colnames(data_NASA_Effort_ReTest) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")



NASA_Effort_long <- pivot_longer(NASA_Effort, cols=c("aclow","bcchall","cchigh","dmlow","emchall","fmhigh"))
div <- rep(c("Mental","Motor"), each = 3)
NASA_Effort_long$div <- rep(div,nrow(NASA_Effort)) 

NASA_Effort_long = split(NASA_Effort_long, NASA_Effort_long$div)

NASA_Effort_coglong = NASA_Effort_long$Mental[,1:3]

NASA_Effort_motlong =  NASA_Effort_long$Motor[,1:3]

###

plot <- ggplot(NASA_Effort_long, aes(name, value, ID, fill = name)) +
  geom_boxplot() +
  labs(x = 'Intensity Level',
       title = 'NASA_Effort',
       y = "NASA_Effort") +
  theme(plot.title = element_text(hjust = 0.5, size = 20),
        legend.title=element_blank(),
        panel.background = element_blank(),
        panel.grid.major.y = element_line(colour = "grey", size = 0.5),
        axis.title = element_text(size = 13),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        strip.text.x = element_text(size = 13)) +
  scale_x_discrete(label = element_blank()) +
  scale_y_continuous(limits = c(-0.015,0.223)) +
  guides(fill="none") +
  scale_fill_manual(values=c("steelblue4", "lightblue1", "cadetblue3", "steelblue4", "lightblue1", "cadetblue3")) +
  annotate(geom = "text", x = 0.83, y = -0.009, label = "Very Low", hjust = 0, vjust = 0, size = 3, color = "gray50") +  
  annotate(geom = "text", x = 0.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.009, label = "Challenging", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.85, y = -0.009, label = "Very High", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.87, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  facet_grid(. ~ div, scales = "free", space = "free")


#ggsave("HRVplot.png", plot, width = 10.5, height = 6.92)

###Statistics

###Test Normality
Norm.NASA_Effort.CogEasy = shapiro.test(NASA_Effort$aclow)
Norm.NASA_Effort.CogMiddle = shapiro.test(NASA_Effort$bcchall)
Norm.NASA_Effort.CogDifficult = shapiro.test(NASA_Effort$cchigh)

Norm.NASA_Effort.MotEasy = shapiro.test(NASA_Effort$dmlow)
Norm.NASA_Effort.MotMiddle = shapiro.test(NASA_Effort$emchall)
Norm.NASA_Effort.MotDifficult = shapiro.test(NASA_Effort$fmhigh)

###Perform Test and Post-Hoc
### Cognitive 
if (Norm.NASA_Effort.CogEasy$p.value > 0.05 & Norm.NASA_Effort.CogMiddle$p.value > 0.05 & Norm.NASA_Effort.CogDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Effort_coglong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Effort_coglong$value, NASA_Effort_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Effort_coglong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Effort_coglong$value, NASA_Effort_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}

###Motor
if (Norm.NASA_Effort.MotEasy$p.value > 0.05 & Norm.NASA_Effort.MotMiddle$p.value > 0.05 & Norm.NASA_Effort.MotDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Effort_motlong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Effort_motlong$value, NASA_Effort_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Effort_motlong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Effort_motlong$value, NASA_Effort_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}


#ICC analysis

## NASA_Effort

## Cognitive

### Easy
data_NASA_Effort_Cog_easy = data.table(data_NASA_Effort_Test$NASA_CognitionEasy, data_NASA_Effort_ReTest$NASA_CognitionEasy)
#data_NASA_Effort_easy_RRn = data_NASA_Effort_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Effort_Cog_easy =ICC(data_NASA_Effort_Cog_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Effort_Cog_easy<- bland.altman.plot(data_NASA_Effort_Test$NASA_CognitionEasy,data_NASA_Effort_ReTest$NASA_CognitionEasy, main="Test-Retest NASA_Effort Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Effort_Cog_easy)

###Middle
data_NASA_Effort_Cog_middle = data.table(data_NASA_Effort_Test$NASA_CognitionOptimal, data_NASA_Effort_ReTest$NASA_CognitionOptimal)
#data_NASA_Effort_Cog_middle = data_NASA_Effort_Cog_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Effort_Cog_middle =ICC(data_NASA_Effort_Cog_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Effort_Cog_middle <-bland.altman.plot(data_NASA_Effort_Test$NASA_CognitionOptimal,data_NASA_Effort_ReTest$NASA_CognitionOptimal, main="Test-Retest NASA_Effort Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Effort_Cog_middle)


###Difficult
data_NASA_Effort_Cog_difficult = data.table(data_NASA_Effort_Test$NASA_CognitionDifficult, data_NASA_Effort_ReTest$NASA_CognitionDifficult)
#data_NASA_Effort_Cog_difficult = data_NASA_Effort_Cog_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Effort_Cog_difficult =ICC(data_NASA_Effort_Cog_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Effort_Cog_difficult<- bland.altman.plot(data_NASA_Effort_Test$NASA_CognitionDifficult,data_NASA_Effort_ReTest$NASA_CognitionDifficult, main="Test-Retest NASA_Effort Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Effort_Cog_difficult)

## Motor

### Easy
data_NASA_Effort_Mot_easy = data.table(data_NASA_Effort_Test$NASA_MotorEasy, data_NASA_Effort_ReTest$NASA_MotorEasy)
#data_NASA_Effort_easy_RRn = data_NASA_Effort_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Effort_Mot_easy =ICC(data_NASA_Effort_Mot_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Effort_Mot_easy<- bland.altman.plot(data_NASA_Effort_Test$NASA_MotorEasy,data_NASA_Effort_ReTest$NASA_MotorEasy, main="Test-Retest NASA_Effort Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Effort_Mot_easy)

###Middle
data_NASA_Effort_Mot_middle = data.table(data_NASA_Effort_Test$NASA_MotorOptimal, data_NASA_Effort_ReTest$NASA_MotorOptimal)
#data_NASA_Effort_Mot_middle = data_NASA_Effort_Mot_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Effort_Mot_middle =ICC(data_NASA_Effort_Mot_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Effort_Mot_middle <-bland.altman.plot(data_NASA_Effort_Test$NASA_MotorOptimal,data_NASA_Effort_ReTest$NASA_MotorOptimal, main="Test-Retest NASA_Effort Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Effort_Mot_middle)


###Difficult
data_NASA_Effort_Mot_difficult = data.table(data_NASA_Effort_Test$NASA_MotorDifficult, data_NASA_Effort_ReTest$NASA_MotorDifficult)
#data_NASA_Effort_Mot_difficult = data_NASA_Effort_Mot_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Effort_Mot_difficult =ICC(data_NASA_Effort_Mot_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Effort_Mot_difficult<- bland.altman.plot(data_NASA_Effort_Test$NASA_MotorDifficult,data_NASA_Effort_ReTest$NASA_MotorDifficult, main="Test-Retest NASA_Effort Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Effort_Mot_difficult)




########################################################################################################################################################


###NASA_Frustration

data_NASA_Frustration_Test = bind_cols(data$Record.ID, data$Frustration, data$Frustration.1 ,data$Frustration.2 ,data$Frustration.3 ,data$Frustration.4 , data$Frustration.5)

data_NASA_Frustration_ReTest = bind_cols(data$Record.ID, data$Frustration.6, data$Frustration.7 ,data$Frustration.8 ,data$Frustration.9 ,data$Frustration.10 , data$Frustration.11)

data_NASA_Frustration_ReTest[10,] = data_NASA_Frustration_Test[10,]

NASA_Frustration = (data_NASA_Frustration_Test+data_NASA_Frustration_ReTest)/2

colnames(NASA_Frustration) <- c("ID","aclow","bcchall","cchigh","dmlow","emchall","fmhigh")
colnames(data_NASA_Frustration_Test) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")
colnames(data_NASA_Frustration_ReTest) <- c("ID","NASA_CognitionEasy","NASA_CognitionOptimal","NASA_CognitionDifficult","NASA_MotorEasy","NASA_MotorOptimal","NASA_MotorDifficult")



NASA_Frustration_long <- pivot_longer(NASA_Frustration, cols=c("aclow","bcchall","cchigh","dmlow","emchall","fmhigh"))
div <- rep(c("Mental","Motor"), each = 3)
NASA_Frustration_long$div <- rep(div,nrow(NASA_Frustration))

NASA_Frustration_long = split(NASA_Frustration_long, NASA_Frustration_long$div)

NASA_Frustration_coglong = NASA_Frustration_long$Mental[,1:3]

NASA_Frustration_motlong =  NASA_Frustration_long$Motor[,1:3]

###

plot <- ggplot(NASA_Frustration_long, aes(name, value, ID, fill = name)) +
  geom_boxplot() +
  labs(x = 'Intensity Level',
       title = 'NASA_Frustration',
       y = "NASA_Frustration") +
  theme(plot.title = element_text(hjust = 0.5, size = 20),
        legend.title=element_blank(),
        panel.background = element_blank(),
        panel.grid.major.y = element_line(colour = "grey", size = 0.5),
        axis.title = element_text(size = 13),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        strip.text.x = element_text(size = 13)) +
  scale_x_discrete(label = element_blank()) +
  scale_y_continuous(limits = c(-0.015,0.223)) +
  guides(fill="none") +
  scale_fill_manual(values=c("steelblue4", "lightblue1", "cadetblue3", "steelblue4", "lightblue1", "cadetblue3")) +
  annotate(geom = "text", x = 0.83, y = -0.009, label = "Very Low", hjust = 0, vjust = 0, size = 3, color = "gray50") +  
  annotate(geom = "text", x = 0.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.009, label = "Challenging", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 1.85, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.85, y = -0.009, label = "Very High", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  annotate(geom = "text", x = 2.87, y = -0.015, label = "Intensity", hjust = 0, vjust = 0, size = 3, color = "gray50") +
  facet_grid(. ~ div, scales = "free", space = "free")


#ggsave("HRVplot.png", plot, width = 10.5, height = 6.92)

###Statistics

###Test Normality
Norm.NASA_Frustration.CogEasy = shapiro.test(NASA_Frustration$aclow)
Norm.NASA_Frustration.CogMiddle = shapiro.test(NASA_Frustration$bcchall)
Norm.NASA_Frustration.CogDifficult = shapiro.test(NASA_Frustration$cchigh)

Norm.NASA_Frustration.MotEasy = shapiro.test(NASA_Frustration$dmlow)
Norm.NASA_Frustration.MotMiddle = shapiro.test(NASA_Frustration$emchall)
Norm.NASA_Frustration.MotDifficult = shapiro.test(NASA_Frustration$fmhigh)

###Perform Test and Post-Hoc
### Cognitive 
if (Norm.NASA_Frustration.CogEasy$p.value > 0.05 & Norm.NASA_Frustration.CogMiddle$p.value > 0.05 & Norm.NASA_Frustration.CogDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Frustration_coglong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Frustration_coglong$value, NASA_Frustration_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Frustration_coglong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Frustration_coglong$value, NASA_Frustration_coglong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}

###Motor
### Cognitive 
if (Norm.NASA_Frustration.MotEasy$p.value > 0.05 & Norm.NASA_Frustration.MotMiddle$p.value > 0.05 & Norm.NASA_Frustration.MotDifficult$p.value > 0.05){
  print(anova_test(data=NASA_Frustration_motlong, dv = value, wid=ID, within = name))
  pairwise.t.test(NASA_Frustration_motlong$value, NASA_Frustration_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
} else { 
  print(friedman.test(data=NASA_Frustration_motlong, value ~ name|ID))
  pairwise.wilcox.test(NASA_Frustration_motlong$value, NASA_Frustration_motlong$name, p.adjust.method = 'bonferroni', paired = TRUE)
}


#ICC analysis

## NASA_Frustration

## Cognitive

### Easy
data_NASA_Frustration_Cog_easy = data.table(data_NASA_Frustration_Test$NASA_CognitionEasy, data_NASA_Frustration_ReTest$NASA_CognitionEasy)
#data_NASA_Frustration_easy_RRn = data_NASA_Frustration_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Frustration_Cog_easy =ICC(data_NASA_Frustration_Cog_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Frustration_Cog_easy<- bland.altman.plot(data_NASA_Frustration_Test$NASA_CognitionEasy,data_NASA_Frustration_ReTest$NASA_CognitionEasy, main="Test-Retest NASA_Frustration Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Frustration_Cog_easy)

###Middle
data_NASA_Frustration_Cog_middle = data.table(data_NASA_Frustration_Test$NASA_CognitionOptimal, data_NASA_Frustration_ReTest$NASA_CognitionOptimal)
#data_NASA_Frustration_Cog_middle = data_NASA_Frustration_Cog_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Frustration_Cog_middle =ICC(data_NASA_Frustration_Cog_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Frustration_Cog_middle <-bland.altman.plot(data_NASA_Frustration_Test$NASA_CognitionOptimal,data_NASA_Frustration_ReTest$NASA_CognitionOptimal, main="Test-Retest NASA_Frustration Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Frustration_Cog_middle)


###Difficult
data_NASA_Frustration_Cog_difficult = data.table(data_NASA_Frustration_Test$NASA_CognitionDifficult, data_NASA_Frustration_ReTest$NASA_CognitionDifficult)
#data_NASA_Frustration_Cog_difficult = data_NASA_Frustration_Cog_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Frustration_Cog_difficult =ICC(data_NASA_Frustration_Cog_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Frustration_Cog_difficult<- bland.altman.plot(data_NASA_Frustration_Test$NASA_CognitionDifficult,data_NASA_Frustration_ReTest$NASA_CognitionDifficult, main="Test-Retest NASA_Frustration Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Frustration_Cog_difficult)

## Motor

### Easy
data_NASA_Frustration_Mot_easy = data.table(data_NASA_Frustration_Test$NASA_MotorEasy, data_NASA_Frustration_ReTest$NASA_MotorEasy)
#data_NASA_Frustration_easy_RRn = data_NASA_Frustration_easy_RRn[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Frustration_Mot_easy =ICC(data_NASA_Frustration_Mot_easy, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Frustration_Mot_easy<- bland.altman.plot(data_NASA_Frustration_Test$NASA_MotorEasy,data_NASA_Frustration_ReTest$NASA_MotorEasy, main="Test-Retest NASA_Frustration Easy", xlab="Means", ylab="differences")

print(ICC_NASA_Frustration_Mot_easy)

###Middle
data_NASA_Frustration_Mot_middle = data.table(data_NASA_Frustration_Test$NASA_MotorOptimal, data_NASA_Frustration_ReTest$NASA_MotorOptimal)
#data_NASA_Frustration_Mot_middle = data_NASA_Frustration_Mot_middle[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Frustration_Mot_middle =ICC(data_NASA_Frustration_Mot_middle, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Frustration_Mot_middle <-bland.altman.plot(data_NASA_Frustration_Test$NASA_MotorOptimal,data_NASA_Frustration_ReTest$NASA_MotorOptimal, main="Test-Retest NASA_Frustration Optimal", xlab="Means", ylab="differences")

print(ICC_NASA_Frustration_Mot_middle)


###Difficult
data_NASA_Frustration_Mot_difficult = data.table(data_NASA_Frustration_Test$NASA_MotorDifficult, data_NASA_Frustration_ReTest$NASA_MotorDifficult)
#data_NASA_Frustration_Mot_difficult = data_NASA_Frustration_Mot_difficult[-3,] #delete unwanted ID (1, 12, 14?)

ICC_NASA_Frustration_Mot_difficult =ICC(data_NASA_Frustration_Mot_difficult, missing=TRUE, alpha = 0.05, lmer = TRUE, check.keys = FALSE) #last 2 variables not needed, they are so by default

BAP_NASA_Frustration_Mot_difficult<- bland.altman.plot(data_NASA_Frustration_Test$NASA_MotorDifficult,data_NASA_Frustration_ReTest$NASA_MotorDifficult, main="Test-Retest NASA_Frustration Difficult", xlab="Means", ylab="differences")

print(ICC_NASA_Frustration_Mot_difficult)

